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Summary 

)  Maximum  entropy  signal  analysis  (MESA)  provides  a  method  of  deriving 
the  auto  power  spectrum  from  a  discrete  number  of  time  series  data;  the 
method  is  reputed  to  be  superior  to  Past  Pourier  Transform  (PPT) 
analysis  When  applied  to  short  data  records  in  that  better  resolution  is 
obtained  and  sidebands  are  reduced. 


In  order  to  investigate  the  efficacy  of  this  method  MESA  has  been 
Implemented  on  computer  controlled  spectrum  analysis  systems  and  a 
number  of  simple  spectra  have  been  examined.  It  is  concluded  that  MESA 
offers  no  advantages  when  large  data  records  are  available  since  the 
quality  of  the  spectra  is  poor  and  long  computation  times  are  required. 
Furthermore  the  method  produces  highly  variable  absolute  levels.  However 
if  only  short  data  records  are  available  then  provided  a  suitable 
sampling  frequency  is  chosen  MESA  does  provide  a  reasonable  means  of 
producing  the  power  spectrum. 
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GLOSSARY  OP  TERMS 


(n+1  )the  coefficient  of  the  m+1  long  prediction  error  filter. 
E  Expectation  value 

f  Frequency 

FPE  Pinal  Prediction  Error 

H  Entropy 

N  Number  of  data  points 

M+l  Order  of  prediction  error  filter 

P(f)  Power  spectrum 

Pm  Power  output  by  (ari-l)th  order  prediction  error  filter 

R^  Measured  values  of  auto-correlation  function 

R(t)  True  values  of  auto-correlation  function 

S(  f )  Linear  spectrum 

t  Tine 

T  Tine  period  of  data  sample 

At  Tine  period  between  data  points 

U(t)  weighting  function 

x<t )j(xi)  Time  series  data  (digitised) 
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1.  INTRODUCTION 


Spectral  data  is  generally  analysed  by  means  of  the  auto-power 
spectrum  P(f)  of  the  signal.  This  is  conventionally  determined  from  the 
measured  time  signal  x(t)  by  means  of  the  Fourier  Transform! 


TIME 

SPECTRUM 

X(t) 


AUTO  CORRELATION  FUNCTION 


lim  l  r 

R<t>  ’ti,  i  ]  x(T)x(T+t)  dr 


i)  PILTER  W(t) 

ii)  FOURIER  TRANSFORM 


i)  FILTER  W(t) 

ii)  FOURIER  TRANSFORM 


LINEAR  SPECTRUM  — 
S(f)  -  F[x(t)] 

-  F[X(t)W(t)] 
in  practice. 


AUTO  POWER  SPECTRUM 
P(f)  -  F[R(  t  )] 

-  P[R(  t  )W(  t )  ]  in 

-  S(f)s*(f) 


practice 


In  practice  the  time  series  data  are  measured  within  some  finite 
time  interval  and  the  analysis  assumes  data  outside  this  interval  to  be 
zero.  To  reduce  the  error  introduced  by  this  assumption  the  available 
data  is  generally  treated  by  some  weighting  function  W(t)  such  as  the 
Hanning  function. 

Maximum  entropy  spectral  analysis  (MESA)  provides  an  alternative 
method  for  determining  the  power  spectrum  P( f )  from  the  available  data 
without  such  modification.  In  principle  the  method  provides  a  means  of 
extending  the  sampled  data  x^  ( iffl )  beyond  the  measurement  window  ( i . e . 
to  i>N)  by  taking  a  weighted  sum  of  the  known  data.  The  weighting 
function  is  known  as  the  Prediction  Error  Filter  and  (as  shall  be 
explained  fully  later)  is  determined  by  minimising  the  error  power 
output  from  a  least  mean  square  fit  of  the  measured  and  'predicted' 
data. 


The  method  was  first  developed  by  Burg  [1]  and  involves  the 
maximisation  of  the  'entropy'  H  of  the  system  with  respect  to  the  power 
spectrum  subject  to  constraints  Imposed  by  the  measured  data.  The  method 
is  reported  to  have  had  much  success  in  improving  the  power  spectrum  in 
that  i)  better  resolution  is  obtained  ii)  sidebands  are  reduced  and  iii) 
more  realistic  power  estimates  axe  produced,  particularly  for  short  data 
records  ( see  for  example  [ l ]  and  [ 2 ] ) . 

In  order  to  investigate  the  claims  and  assess  ths  potential  of  the 
method  for  use  on  machinery  noise  data  MESA  has  been  implemented  on  two 
analysis  systems  (the  GENRAD  2508  Signal  analysis  system  and  the  Hewlett 
Packard  5423A  Spectral  dynamics  analyser  used  un  conjunction  with  the  HP 
9845  deftk  top  computer)  and  has  been  used  to  analyse  several  simple 
spectra. 


2.  THEORY 


The  concept  of  entropy  as  used  in  MESA  was  originally  derived  from 
the  field  of  information  theory.  The  entropy  of  a  system  is  the  average 
rate  of  information  conveyed  and  is  given  by 


n 

H  -  E  P,  ln(P  )  1 

i-t 

where  p^  is  the  probability  that  the  ith  event  of  n  possible  events 
occurs . 


Prom  equation  1  it  can  be  seen  that  the  entropy  is  minimised  when 
the  probabilities  are  equal,  i.e.  when  the  system  is  most  ordered. 

The  application  of  the  concept  of  entropy  to  spectral  analysis  therefore 
involves  the  minimisation  of  the  entropy  for  a  stationary  Gaussian 
process  which  has  been  shown  [1,2,3]  to  be  proportional  to 


H  -  J  In  P(f)  df  2 

Where  P( f )  is  the  power  spectrum  and  must  be  consistent  with  the 
measured  values  of  the  auto-correlation  function  ( assuming  these  to  be 
known  exactly). 

The  minimisation  procedure  is  described  in  Appendix  1  and  results  in 
the  following  expression  for  the  power  spectrum; 


P(f)  - 


m 


At 


exp(  -2rrifnAt ) 


2 
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where  At  is  the  time  interval  between  sampled  data  points,  Pm  is  the 
power  output  by  the  prediction  error  filter  of  order  hh-1  having  elements 
a^,,;  these  elements  are  the  solutions  to  the  equation 
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where  R*  are  the  measured  values  of  the  auto-correlation  function. 


This  analysis  is  performed  on  the  assumption  that  the  values 
Rfc  (  -  ri--)  “  ”  E[xi  Xj]  )  of  the  auto  correlation  function  axe 

exactly  known.  However,  in  practice  R*  are  determined  from  truncated 
sampled  time  series  data  x^  whence  Ry  are  no  longer  exact  but  are  given 
by 


M-k 


E 

i-l 


xi  Vk 
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Mow,  by  considering  the  final  row  of  4  it  can  be  seen  that  Ry  can  be 
expressed  in  terms  of  the  prediction  error  filter  coefficients  and 
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V 


proceeding  values  of  the  auto-correlation  coefficients: 


\  - 


Viai 


Burg  [4]  proposed  a  means  of  obtaining  more  realistic  values  of  R* 
by  iterating  values  for  the  filter  coefficients  by  successive 
incrementation  of  the  dimension  of  the  filter.  (It  is  convenient  at 
this  stage  to  redefine  the  elements  a-j  as  a^j  where  m  refers  to  the 
number  of  iterations  performed . )  This  method  is  based  essentially  on  a 
least  squares  fitting  of  the  true  data  to  the  output  of  the  prediction 
error  filter.  Consider,  for  example,  the  two  point  prediction  error 
filter 


In  effect,  this  filter  serves  to  •predict'  the  (i+l)th  value  of  the 
data  xi+1( pred )  from  the  ith  measured  value  x^  by  setting  xi+1(pred)  » 
allxi<  In  order  to  obtain  the  ’best*  value  of  alt  the  filter  is  scanned 
along  the  data  sample  and  a  least  squares  fit  is  performed  between  the 
measured  and  predicted  values  to  determine  the  value  of  an  which 
minimises  the  'error*  P1.  In  order  to  ensure  that  the  resulting  filter 
is  stable,  it  is  required  that  the  coefficients  a^  should  be  less  than 
unity.  To  satisfy  this  criterion  Burg  [4]  suggested  that  the  forward 
power  output  of  the  filter  should  be  averaged  with  the  power  output  by 
the  filter  operated  in  reverse: 
reverse. 


Forward  filter 

au  1 

x. . .pred  -  a. _x. 

au  1 

1  an  < 


Backward  filter 
xiPred  -  axlxi+1 


1  a. 


1  *  r  2  2  1 

-i ts=T)  (  xi+r  auxi)  +  (xrauxi+i>  }  7 


Minimising  P„  with  respect  to  a. . :  - — 

x  xi  da , 


whence 


ail  " 


2  lmt  xixitl 


M— 1 

U  *i  +  Xi+1> 
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the  prediction  error  filter  could  become  more  accurate  on  increasing 
the  number  of  data  points  incorporated  in  the  estimate  resulting  in  am 
increased  filter  length.  The  procedure  for  calculating  a  new  set  of 
filter  coefficients  on  incrementing  the  filter  length  is  described  in 
Appendix  2.  In  general  we  find  the  coefficients  through  the  recursion 
relation 


ode 


a  .  .  -  a 
m-l  k  a 


m-l  m-k 


k<m 


with  a^,,  being  derived  from  proceeding  coefficients  and  the  original 
data 

A  recursion  for  the  error  power  output  can  also  be  derived  (Appendix 
3)  and  is  given  by 


P  -  P  ,  (  1  -  a 
m  m-l  ■ 

where,  from  4  and  S 
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P_  -  R„  - 
0  0 


N 


i-l 
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Returning  to  the  expression  for  the  power  spectrum,  equation  3, it 
can  be  seen  that  knowing  the  elements  a^-j  and  the  error  power  output  P^ 
the  explicit  derivation  of  the  auto-correlation  coefficients  becomes 
unnecessary  and  solution  of  4  and  6  can  be  by-passed. 


*  tw  wntmc non  intnnn  ptt/pkw 

The  choice  of  length  for  the  prediction  error  filter  still  appears 
to  be  a  matter  of  concern.  Too  short  a  filter  results  in  smoothing  of 
the  data  while  too  long  a  filter  causes  unwanted  splitting  of  the  peaks, 
three  criteria  are  available  for  choosing  the  length  of  the  filter 
[6,7]: 

DflnajL  Piadlgtion  Error  (PPE)  Criterion 

* 

This  procedure  was  proposed  by  Akaiks  [8,9]  and  consists  of 
determining  the  filter  length  M  for  which  the  FPE  -  E((x^-x^)2}  is 
minimised.  It  is  shewn  that 

(N+m+l ) 

PPE(m)  -  -  P  12 

(N-m-1)  “ 

should  be  minimised. 

ii)  Infpxroaupn  JThcprttlc  £iltexign  (AIK) 

this  is  based  on  the  minimisation  of  the  log  liklihood  of  the 
prediction  error  variance 


iii) 


This  ia  based  on  minimising  the  difference  between  the  mean  square 
errors  of  the  estimated  and  true  filters 
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The  first  of  these  methods,  viz  the  PPE  criterion,  has  been  most 
widely  used  in  the  literature  and  was  adopted  in  this  assessment. 


Olrych  and  Bishop  [2]  recommend  M-N/2  as  a  generally  suitable  filter 
length;  however,  from  the  point  of  view  of  computer  time  a  shorter 
filter  is  preferable. 


4  COMPUTER  PROGRAMS 

MESA  has  been  implemented  on  both  the  GENRAD  2508  Signal  analysis 
system  (Time  Series  Language,  TSL)  and  on  the  HP9845  desk  top  computer 
(BASIC)  used  in  conjunction  with  the  HP  5423A  Spectral  Dynamics 
Analyser.  Figure  1  shows  the  overall  structure  of  both  programs  while 
Figure  2  shows  the  algorithm  used  for  determining  the  power  output  Pm 
and  the  filter  coefficients  In  both  causes  the  FPE  is  plotted  as  a 

function  of  m  to  allow  the  user  to  assess  the  length  of  the  filter  used; 
this  length  can  be  changed  if  required  before  going  on  to  determine  the 
power  spectrum. 

On  running  the  program  on  the  GENRAD  the  levels  of  the  power  spectra 
exhibited  high  variability;  this  was  initially  attributed  to  rounding 
errors  incorporated  in  the  TSL  software.  Since  the  method  is 
particularly  sensitive  to  such  errors  the  majority  of  plots  presented  in 
this  study  were  derived  using  the  Hewlett  Packard  system  where  the 
increased  precision  available  would  minimise  their  effect.  However  (as 
described  in  the  next  section)  the  amplitudes  were  still  found  to  vary 
considerably. 


5.  AMMiWXa  PERFORMANCE 

The  MESA  method  was  tested  on  various  signals  including 
i )  Pure  sine  wave  ( lv  rms ,  50Rz ) 

11)  Two  adjacent  sine  waves  (lv  rms,  50Kz;  lv  rms,  49Rz) 
ill )  Square  wave  ( lv  rms ,  lOKz ) 
iv)  Typical  machinery  noise  . 

Ply  SOB*  am m.  wart 

This  simple  case  was  explored  fairly  extensively  with  the  sampling 
frequency,  number  of  data  points  and  number  of  iterations  being  varied. 

The  general  observations  from  this  investigation  were  as  follows; 

a)  A  two  point  error  filter  (M»l)  was  incapable  of  producing  a 
suitable  spectrum. 
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to)  One  iteration  ( M-2 )  was  sufficient  to  produce  a  peak  in  the 
spectrum;  increasing  the  number  of  iterations  did  not  affect  the 
position  of  this  peak  but  did  improve  the  resolution  and  dynamic  range. 
This  conclusion  is  illustrated  in  Figure  3  which  shows  spectra  derived 
from  2,8  and  32  iterations  on  64  data  points  sampled  at  256Hz.  Note  that 
for  M"32(*N/2)  peak  splitting  occurred;  such  splitting  for  long  filters 
has  been  reported  in  the  literature  (10].  The  corresponding  variation  of 
FPE  with  M  is  shown  in  Figure  4:  there  is  no  obvious  minimum  but  a 
levelling  off  occurs  at  around  M=8  and  this  does  correspond  to  a 
reasonable  spectrum. 

c)  The  position  and  width  of  the  peak  varied  with  the  number  of  data 
points  and  with  the  sampling  frequency.  This  is  illustrated  in  Figure  5 
which  shows  the  peak  position  converging  on  50Hz  as  the  number  of  points 
in  the  sample  increased  (sampling  frequency  3200Hz,  M=2 ) .  on  decreasing 
the  sampling  frequency  the  convergence  became  more  rapid  ( Figure  6 )  such 
that  a  true  spectrum  could  be  produced  with  fewer  points.  Indeed,  at  the 
Nyquist  frequency  an  accurate  (excepting  peak  height)  spectrum  could  be 
produced  with  as  few  as  three  data  points  ( Figure  7 ) ;  However  the 
reflection  of  the  zero  frequency  peak  at  100Hz  shows  that  aliasing  was 
occurring . 

The  resolution  also  improved  as  the  peak  converged  to  the  correct 
frequency. 

d)  The  method  was  incapable  of  producing  the  correct  absolute  level  of 
the  spectrum. 

The  major  conclusions  from  this  simple  investigation  were  that  mesa 
cannot  be  used  to  give  accurate  levels.  Nethertheless  it  is  capable  of 
producing  a  well  resolved  spectrum  (of  a  sine  wave)  using  very  few  data 
points  providing  the  sampling  frequency  is  low;  however  on  approaching 
the  tfyquist  frequency  aliasing  is  seen  to  occur. 

ii>  Mjacetrt .  9i,n<t  sasaa-1 

Figure  8  shows  the  ME  spectra  derived  from  a  128  point  data  set 
sampled  at  200Hz  to  avoid  aliasing.  It  can  be  seen  that  four  iterations 
were  insufficient  to  resolve  the  peaks  while  8  and  32  iterations 
resolved  two  (unequal)  peaks  at  around  49Hz  and  51Hz,  the  latter 
frequency  being  incorrect.  On  decreasing  the  length  of  the  data  sample 
to  16  points  MESA  was  still  capable  of  resolving  two  peaks  (albeit  at 
46Hz  and  SIHz)  when  8  iterations  were  used;  however  an  8  point  data 
sample  was  insufficient  for  these  to  be  resolved. 

The  FPE  criterion  produced  similar  results  as  for  the  single  sine 
wave  and  did  not  provide  any  guidance  as  to  the  best  filter  length. 


in  wa.v9 

In  order  to  assess  MESA'S  ability  to  cope  with  harmonics,  a  10Hz 
square  wave  was  examined;  the  majority  of  tests  were  performed  with  a 
sampling  frequency  of  200 Hz  in  order  to  examine  harmonics  up  to  100Hz. 

Figure  9  shows  the  analysis  results  for  a  128  point  data  sample 
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using  different  filter  lengths.  Eight  iterations  were  sufficient  to 
produce  a  fair  representation  of  the  spectrum  in  that  the  first  four  of 
the  odd  harmonics  occurred  close  to  the  correct  frequencies  although  the 
peak  heights  were  subject  to  variability.  The  fifth  odd  harmonic  was 
only  resolved  in  the  case  of  the  long  filter,  increasing  the  number  of 
iterations  improved  the  resolution.  For  m-32  peaks  were  observed  at  even 
harmonics ;  since  these  were  also  observed  on  performing  a  PFT  analysis 
they  were  attributed  to  to  the  imperfection  of  the  input  sine  wave. 

Application  of  the  PPE  criterion  proved  unhelpful  in  choosing  the 
•best’  filter  length  to  be  used;  after  32  iterations  no  true  minimum  was 
attained  (Pigure  10).  (Longer  filters  were  not  considered  due  to  the 
impracticality  of  computation  time  -  more  them  ten  minutes  were  required 
to  produce  a  spectrum  using  N*=128,  M =32)  Pigure  11  shows  the  ME  spectra 
of  shorter  length  samples  ( M=8  in  each  case)  and  shows  that  even  with 
only  16  data  points  a  representative  spectrum  can  be  produced.  However 
Pigure  12  indicates  the  importance  of  choosing  a  suitable  sampling 
frequency:  Too  low  a  sampling  frequency  caused  aliasing  (fa=l00Rz) 
while  at  too  high  a  frequency  the  small  data  sample  only  represented  a 
fraction  (<l/2)  of  a  cycle. 

The  conclusions  drawn  from  this  test  are  that 

i )  MESA  is  capable  of  producing  a  representative  spectrum  ( in  terms  of 
peak  position)  of  the  first  few  harmonics  of  a  10Hz  square  wave  even 
when  using  a  small  data  sample.  However  higher  harmonics  were  not 
resolved. 

it)  The  performance  is  critically  sensitive  to  the  sampling  frequency. 

iii)  The  resolution  can  be  improved  by  using  larger  data  samples  but 
computation  time  becomes  impracticable. 

iv)  Peak  levels  are  incorrect  even  in  relative  terms. 

v)  The  PPE  criterion  offers  little  guidance  in  choosing  filter  length. 


iv)  Machinery. 

A  30kW  motor  was  chosen  to  represent  a  simple 
example  of  machinery  noise;  the  'true'  spectrum  ( FFT  using  512  data 
points )  for  a  random  run  is  shown  in  Pigure  13  ( not  corrected  for 
measurement  parameters). 

The  MESA  performance  on  the  first  128  of  these  data  points  is  shown 
in  Pigure  14  for  64  (N/2  as  recommended  in  [2]),  20  (where  PPE  levelled 
off),  and  5  iterations  of  the  prediction  error  filter.  It  can  be  seen 
that  with  only  5  iterations  MESA  was  able  to  identify  the  two  major 
peaks  while  with  M-20  and  M-64  minor  peaks  were  also  resolved.  However 
the  performance  above  about  70  Hz  was  disappointing,  a  feature  also 
observed  on  examination  of  the  square  wave.  Pigure  16  illustrates  the 
performance  of  MESA  using  only  16  data  points  and  shows  that  the  major 
peaks  could  be  accurately  resolved  although  minor  peaks  had  disappeared 
due  to  lack  of  information.  A  further  test  using  only  8  data  points 
succeeded  in  resolving  the  major  ( 30Hz )  peak  accurately  and  exhibited  a 
second  peak  at  around  70Hz  (Instead  of  60Hz) 


6 


A  short  program  was  Implemented  on  the  GENRAD  2508  to  determine  the 
auto-power  spectrum  by  applying  the  FFT  to  small  samples  of  data  points. 
In  order  to  preserve  the  resolution  of  the  output  spectrum  zero  packing 
was  utilised  (viz  the  FFT  was  performed  on  a  128  point  array 
incorporating  N  data  points,  all  other  elements  being  set  to  zero). 

Results  for  the  50Hz  sine  wave  and  10Hz  square  wave  are  shown  in 
Figure  17  for  16,  32  and  64  data  points,  both  waves  being  sampled  at 
200HZ . 

Comparing  these  figures  with  those  derived  using  MESA  (Pigures  S,  7, 
and  9)  it  can  be  seen  that  the  MESA  performance  is  by  far  the  more 
satisfactory  when  only  16  or  32  data  points  are  present  in  that  peaks 
can  be  easily  resolved;  for  64  points  MESA  improves  on  the  dynamic  range 
and  resolution. 


7  SIGNAL  VARIANCE  -  PISCffSSIQW 

Referring  to  equation  3  the  ME  power  spectrum  can  be  seen  to  be  the 
result  of 


Error  power  spectrum  level 

p(f)  a  - -  15 

I  F.T.  of  prediction  error  filter  I 

Thus  in  order  for  peaks  to  be  observed  in  P(f)  the  Fourier  transform 
must  exhibit  minima.  The  peak  amplitude  is  therefore  a  result  of  a 
limiting  process  in  which  both  the  error  power  and  Fourier  transform 
tend  to  zero  but  in  such  a  way  that  the  power  spectrum  level  remains 
finite . 

This  probably  explains  in  part  the  large  variation  observed  in  the 
amplitudes.  This  is  also  commented  on  in  [2]  as  being  associated  with 
harmonic  signals. 


A  second  source  of  variance  arises  from  rounding  errors.  These  are 
particularly  troublesome  in  the  determination  of  the  error  power  Pm 
(equation  10)  when  tends  to  unity.  This  would  occur  (see  equation 
A2.6)  when  b£t  which  is  effectively  saying  that  forward  and 

backward  predictors  give  the  same  output  as  would  occur  if  the  sampling 
points  were  symmetric  in  the  signal.  (This  would  only  be  true  of  a 
noiseless  signal  Which  would  be  unlikely  in  practice. ) 


8  C0HCW5I0M8 

MESA  has  been  shown  to  be  a  viable  means  of  analysing  data  provided 
suitable  sampling  frequencies  and  filter  lengths  are  chosen.  However  a 
major  shortcoming  is  that  the  lsvsls  derived  are  extremely  variable. 
Frequency  shifting  can  occur  but  can  be  reduced  by  Choosing  a  sampling 
frequency  close  to  the  Nyqulst  frequency.  Spectral  resolution  is 
improved  by  increasing  the  length  of  the  prediction  error  filter; 
however  if  this  filter  is  too  long  spurious  peaks  can  be  produced. 


j _ I 
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Akaike’s  FPE  criterion  did  not  prove  particularly  helpful  as  a  guide  to 
filter  length  in  any  of  the  examples  considered  -  no  clear  minima 
occured  while  the  value  of  M  following  a  rapid  decrease  in  FPE  seemed  as 
good  a  choice  as  any. 

When  long  data  samples  sure  available  MESA  offers  no  advantages  over 
the  conventional  FFT  means  of  analysis :  the  spectra  sure  of  poor  quality, 
do  not  provide  reliable  absolute  levels,  and  their  execution  is 
unacceptably  slow.  However,  when  only  few  data  points  are  available  MESA 
has  definate  advantages  in  that  major  peaks  can  be  clearly  resolved 
although  reservation  is  reconraended  with  regard  their  position. 

In  the  context  of  time  series  analysis  of  machinery  noise  MESA  is 
unlikely  to  be  of  practicable  value  since  data  samples  can  generally  be 
made  sufficiently  long  to  allow  conventional  analysis.  However  for 
problems  involving  spatial  analysis  where  long  data  samples  are 
inconvenient  if  not  impossible,  MESA  provides  a  potentially  useful  tool. 


L.V.  EMBLING.  (HSO) 


The  entropy  for  a  stationary  Gaussian  process  is  proportional  to 


-  J  ln[p(f)] 


where  the  power  spectrum  P(f)  must  be  consistent  with  the  measured 
values  Rfc  (X=o  to  o)  of  the  auto-correlation  function 


*  |  df  P(f)  exp<  2rrifXat  ] 


with  Rjj  «*  Ri-j  “  Rij  =  E[x^Xj ];  j  =  i+X 

The  procedure  for  minimising  H  with  respect  to  P(f)  in  order  to 
determine  P(f)  is  described  fully  in  [3].  The  method  consists  in 
introducing  Lagrange  multipliers  ©^  where  ©)<=  ©I*  and  demanding  that 

6(  L+H )  «  0  A1.2 


Where 


L-=  ^  s 

JC*— SI 


Combining  A1.1,A1.2  and  AJL.3, 


n 

O  -  S(LtH)  -  Jdf  6(P(f))[  -  £  0*  exp( 2>ri£XAt ) 


and  since  P(f)  is  arbitrary s- 


P(«)  -  [  £  ©£  exp(  2rri£XAt )  ]  1 


Imposing  the  condition  that  P<  f )  is  positive  and  integrable  it  can  be 

expressed  as 

P(f)  “  - - -  KL. 5 

h(f)  A  («) 


where 


A<f)  -  £  yk  exp(  2rrifXAt ) 

X-0 


'  =  yX  yX+n 
X-0 


Multiplying  both  sides  of  A1.5  by  A*exp( 2nifXAt )  and  integrating  in 
the  upper  half  plane  it  is  eventually  shown  by  considering  the 
constraint  equations  Al.l  and  using  the  periodicity  of  A(f)  that 


A.!7!  ”  y„  6X,0 
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or, in  matrix  notation 


[Rlty]  -  l/y0C8] 


A1 .  9 


where 
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fro m  AI .5  and  At. 6 

P(f)  “ 

A<f)l~2  -  |A*(f)l“2 

m-i 

-2 

- 

£  yv  exp(  -2>rifJtAt ) 

Jc-0 

m 

i  [y]T  [«<«)]  i-2 

-  yj  I  <C*rl£*l>*  [«<*)]  I"2 


where 


C«<f)]  - 

Pro*  Al.ll  we  find 


exp<  -2iri(  n-1  )tAt ) 


1  r„,-l 
y0  *  y0CR1U 


whence 


P(f)  - 


P«At 


l<  [Rf1!*!  >T  te< f )]  I2 
Rewriting  AI. 15  in  the  form  ueed  by  [1,5  and  6] 

P_  At 


P<«)  * 


1  -  C  a  exp< -2irif  nAt ) 

n-l  n 


AI .  10 


AI .  11 


A1.12 


AI .  13 


AI.  14 


AI .  15 


1  A 


L  _ 


where  PB  At  -  yo2  And  an  =  -yn/yo  Are  solutions  of  4  (vix  Al.9);  Pm 
is  known  as  the  power  output  of  the  prediction  error  filter,  [an]  and 
[R]  is  the  TOEPLITZ  matrix  [1]: 


P(f)  is  the  MAXIMUM  ENTROPY  solution  to  the  power  spectrum. 


&EPSNPIX  2 

ITERATIVE  PROCEDURE 

The  average  error  power  output  of  A  m+1  long  prediction  error  filter 
C«mlT  ‘[1  «n  •••  i*  given  by  [5]i 


N-m 

r  ®  1 2  | 

III 

| 

E 

t-1  11 

l  V  }  *  1 

.  Xt+m  ^  amkXt-Mn-k  |  j 
k*»i 

A2.1 

The  filter  coefficients  Are  given  by  the  recurrance  relations  (see 
Appendix  3) 


a«k  " 


—  a  a 

-l  k  mm  m-1  or-k 


k«m 


Where  a 


0  for  k>a. 


In  order  to  determine  the  coefficient  a^.  we  insert  10  in  A2.1  to  give 
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A2.3 


Note 


b  ax  and  b  ■  x 
It  t  It  t+l 


A2.4 


Since  b^  and  b^  are  independent  of  ^  the  minimisation  of  Pm 
with  respect  to  a^  gives 


1  •  *2  2 
— - —  r  -4b  b  +  2b  a  +  2b  a 

2(N-m)  “  ,  ntnt  at  mn  mtnm 

and  hence 


-  0 


[MB  this  reduces  to  equation  9  as  required  on  substituting  m»l] 
Furthermore ,  s ince 


d2P 


1  M"  *  •  2  2 
(NHn)  bmt  +  b»t  >  >  0 


A2 .6 


da 

the  extreeui  is  a  minimum  as  required. 

Combining  A2.1  and  A2 . 2  gives  recursion  relations  for  b^  and 


b  -  b  -  a  b 

art  m-1  t  m-1  m-1  m-1  t 


b  ■  b  —  a  b 

at  m-1  t+l  m-l  m-l  m-l  t+l 


A2.7 


ftgEEMD.IX._3 

REOJRSIOM  RELATIONSHIPS 

Consider  the  equation  defining  the  power  output  Pm_1  using  an  m 
point  prediction  error  filter  (cf  eq  4): 
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Increasing  the  order  of  the  matrices  by  1: 
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Where 


R  R„ 

♦  •  R 

R 

1 

p 

o  1 

m-1 

m 

m-l 

R. 

• 

-a 

0 

1  0 

• 

m-l 

1 

•  • 

. 

. 

• 

• 

mt 

• 

R  ,  R 

•  • 

R\ 

-a 

0 

m-1  m 

0 

l 

m-l 

m-l 

R  R 

•  R 

R 

0 

A 

m  m-l 

1 

oj 

m-lj 

m-1 


A  »  R  -  r  r  a 

Vl  m  ^VjVi  j 


A3. 2 


Mddng  use  of  the 
be  reversed  to  give* 


of  A3. 2  the  prediction  error  filter  can 
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A3. 3 


the  solution  for  the  mflth  power  estimate  Pa  is  given  by  equation  4 
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lot  the  prediction  error  filter  be  defined  in  terras  of  the  (m-l)th 
filter: 


1 
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am-l  m-i 
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0 
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. 

Then  equating  elements  on  both  sides  of  A3. 4  we  find 


-a  -  r 


and 


a«*-  Vu  +  r  Vi  «-k 
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m-l  k  m-l  n-k 


The  corresponding  power  output  is  given  by 
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and 
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A3. 5 


Hence 


p— i  -  r2  p, 
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LOG<FPE<H)^FPEU)  POWER  SPECTRUM  dB  ra  l  V#*2 


FREQUENCY  Hz 

FIGURE  3.  VARIATION  OF  MAXIMUM  ENTROPY  SPECTRUM 
WITH  NUMBER  OF  ITERATIONS 
l  iv  rms  50Hz  sins  wave  samp  ted  at  25SHz.) 


8  3  10  IS  20  25  30 

NO  OF  ITERATIONS 

FIGURE  4.  VARIATION  OF  FPE  WITH  M 
(50Hz  Sins  Wavs t  N“64 j  Sampling  Frequency  256Hz) 
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POHER  SPECTRUM  dB  re  1  V»#2  POWER  SPECTRUM  dB  re  1  V**2 


FREQUENCY  Hz 

FIGURE  7  MAXIMUM  ENTROPY  SPECTRUM  OF  50Hz  SINE  WRVE 
SAMPLED  AT  NYQUIST  FREQUENCY 


FREQUENCY  Hz 

FIGURE  8.  MESA  PERFORMANCE  ON  TWO  SINE  WAVES 
lv  raw  st  49Hz  end  50Hz|  Ssnpllng  Frequency  200Hz 
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MB KD  TM30S7 


N-  126 


M-  32 


FREQUENCY  Hz 

FIGURE  9.  MAXIMUM  ENTROPY  SPECTRUM  OF  10Hz  SQUARE  WRVE 
(Sampling  Frequency  208Hz) 
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NUMBER  OF  ITERATIONS 

FIGURE  10.  VARIATION  OF  FPE  WITH  M 

(10Hz  Square  Waves  N»128j  Sampling  Frequency  200Hz) 
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FREQUENCY  Hz 

FIGURE  II.  VRRIRTION  OF  MAXIMUM  ENTROPY  SPECTRUM 
WITH  NUMBER  OF  DRTR  POINTS 
( 10H*  Square  Wave  Sampled  at  200Hz> 


Sampling  Frequency  -  100Hz 

200Hz 

400Hz 


FREQUENCY  Hz 

FIGURE  12.  VRRIRTION  OF  MAXIMUM  ENTROPY  SPECTRUM 
WITH  SRMPLING  FREQUENCY  . 

( 10Hx  Square  Wavei  N*16|  M-0 ) 


FREQUENCY  Hz 
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FIGURE  17  FFT  PERFORMANCE  ON  SMALL  DATA  SAMPLES 


